##########################################
########################################## 
## STUDY 4
##########################################
########################################## 
setwd(getwd())
source("libraries.R")

# load data for study
data_raw <- haven::read_dta(here("Data", "study4.dta"))

##########################################
########################################## 
## Recode Variables
##########################################
########################################## 
df <- data_raw %>% 
  mutate(female = gender,
         race = case_when(
           race == 2 ~ 0,
           TRUE ~ 1
         ),
         
         #Party Identification (1 = Democrat; 0 = Republican)
         repdem = case_when(
           rep1 == 1 | rep1 == 2 | repdem2 == 1 ~ 0,
           dem1 == 1 | dem1 == 2 | repdem2 == 2 ~ 1,
           TRUE ~ NA_real_)) %>% 
  
  #Dummy indicating whether participant is Independent or not
  mutate(independents = case_when(
    repdem2 == 3 ~ 1,
    TRUE ~ NA_real_
  )) %>% 
  
  #Political Violence Intentions
  mutate_at(vars(conflict1_1:conflict2_7), ~na_if(., 6)) %>% 
  mutate_at(vars(conflict1_1:conflict2_7), ~na_if(., 99)) %>% 
  mutate(violentact = rowMeans(select(., c("conflict2_2", "conflict2_3", "conflict2_5", "conflict2_6", "conflict2_7")), 
                               na.rm = TRUE)) %>% 
  
  #Need for Chaos
  mutate(chaos = rowMeans(select(., c("needchaos_1", "needchaos_3", "needchaos_4", 
                                      "needchaos_5", "needchaos_6", "needchaos_9", "needchaos_10")), 
                          na.rm = TRUE)) %>% 
  
  #Political Rumors: Believing and Sharing of DEMocratic and REPublican rumors
  #In main analysis, we focus on reverse coded measures, indicated by an "r"
  mutate(true_DEM = rowMeans(select(., c("dem_true1","dem_true2","dem_true3","dem_true4","dem_true5")), na.rm = TRUE)) %>% 
  mutate(share_DEM = rowMeans(select(., c("dem_share1","dem_share2","dem_share3","dem_share4","dem_share5")), na.rm = TRUE)) %>% 
  mutate(true_REP = rowMeans(select(., c("rep_true1","rep_true2","rep_true3","rep_true4","rep_true5")), na.rm = TRUE)) %>% 
  mutate(share_REP = rowMeans(select(., c("rep_share1","rep_share2","rep_share3","rep_share4","rep_share5")), na.rm = TRUE)) %>% 
  mutate(true_DEMr = rowMeans(select(., c("dem_true_R_1","dem_true_R_2","dem_true_R_3","dem_true_R_4","dem_true_R_5")), na.rm = TRUE)) %>% 
  mutate(share_DEMr = rowMeans(select(., c("dem_share_R_1","dem_share_R_2","dem_share_R_3","dem_share_R_4","dem_share_R_5")), na.rm = TRUE)) %>% 
  mutate(true_REPr = rowMeans(select(., c("rep_true_R_1","rep_true_R_2","rep_true_R_3","rep_true_R_4","rep_true_R_5")), na.rm = TRUE)) %>% 
  mutate(share_REPr = rowMeans(select(., c("rep_share_R_1","rep_share_R_2","rep_share_R_3","rep_share_R_4","rep_share_R_5")), na.rm = TRUE)) %>% 

  #Motivations for sharing hostile rumors about Democrats and Republicans
  mutate(mobilizeDem = rowMeans(select(., c("dem_mob_1", "dem_mob_2", "dem_mob_3" , "dem_mob_4", "dem_mob_5", "dem_mob_6")), 
                                na.rm = TRUE)) %>% 
  mutate(mobilizeRep = rowMeans(select(., c("rep_mob_1", "rep_mob_2", "rep_mob_3", "rep_mob_4", "rep_mob_5", "rep_mob_6")), 
                                na.rm = TRUE)) %>% 
  mutate(trueDem = rowMeans(select(., c("dem_true_1", "dem_true_2", "dem_true_3", "dem_true_4", "dem_true_5", "dem_true_6")), 
                            na.rm = TRUE)) %>% 
  mutate(trueRep = rowMeans(select(., c("rep_true_1", "rep_true_2", "rep_true_3", "rep_true_4", "rep_true_5", "rep_true_6")), 
                            na.rm = TRUE)) %>% 
  mutate(amuseDem = rowMeans(select(., c("dem_amuse_1", "dem_amuse_2", "dem_amuse_3", "dem_amuse_4", "dem_amuse_5", "dem_amuse_6")), 
                             na.rm = TRUE)) %>% 
  mutate(amuseRep = rowMeans(select(., c("rep_amuse_1", "rep_amuse_2", "rep_amuse_3", "rep_amuse_4", "rep_amuse_5", "rep_amuse_6")), 
                             na.rm = TRUE)) %>% 
  mutate(conseqDem = rowMeans(select(., c("dem_conseq_1", "dem_conseq_2", "dem_conseq_3", "dem_conseq_4", "dem_conseq_5", "dem_conseq_6")), 
                              na.rm = TRUE))  %>% 
  mutate(conseqRep = rowMeans(select(., c("rep_conseq_1", "rep_conseq_2", "rep_conseq_3", "rep_conseq_4", "rep_conseq_5", "rep_conseq_6")), 
                              na.rm = TRUE)) %>% 
  
  #Standardize variables
  mutate(across(.cols = c(chaos, violentact,true_DEM, share_DEM, true_REP, share_REP, mobilizeDem, 
                          mobilizeRep, trueDem, trueRep, amuseRep, amuseDem, conseqRep, conseqDem, female, age, race, edu,
                          true_DEMr, share_DEMr, true_REPr, share_REPr),
                ~ (scale(.) %>%  as.vector),
                .names = "{.col}_z_score")) %>% 
  
  #Scale variables to range from 0 to 1
  mutate(across(.cols = c(chaos, violentact,true_DEM, share_DEM, true_REP, share_REP, mobilizeDem, 
                          mobilizeRep, trueDem, trueRep, amuseRep, amuseDem, conseqRep, conseqDem, female, age, race, edu,
                          true_DEMr, share_DEMr, true_REPr, share_REPr),
                ~ (zero1(.) %>% as.vector),
                .names = "{.col}_01")) %>% 
  mutate(mobilizeAll_01 = zero1(rowMeans(select(., c("mobilizeDem_01", "mobilizeRep_01")), na.rm = TRUE)),
         amuseAll_01 = zero1(rowMeans(select(., c("amuseDem_01", "amuseRep_01")), na.rm = TRUE)),
         conseqAll_01 = zero1(rowMeans(select(., c("conseqDem_01", "conseqRep_01")), na.rm = TRUE)),
         trueAll_01 = zero1(rowMeans(select(., c("trueDem_01", "trueRep_01")), na.rm = TRUE)))

##########################################
########################################## 
## Sample Description 
## (See Supplementary Materials SM1, Table S2)
##########################################
########################################## 
#Gender
print(paste(round(mean(df$female_01), 2)*100 , "percent of participants are female"))
#Age
print(paste("Mean age of participants:" , round(mean(df$age) , 2), 
            ", SD:", round(sd(df$age) , 2)))
#Education
print(paste("Less than high school:", round(prop.table(table(df$edu)), 2)[1]*100, "percent"))
print(paste("High school graduate:", round(prop.table(table(df$edu)), 2)[2]*100, "percent"))
print(paste(" Some college, but no degree:", round(prop.table(table(df$edu)), 2)[3]*100, "percent"))
print(paste("2 year college degree:", round(prop.table(table(df$edu)), 2)[4]*100, "percent"))
print(paste("4 year college degree:", round(prop.table(table(df$edu)), 2)[5]*100, "percent"))
print(paste("Graduate or professional degree:", round(prop.table(table(df$edu)), 2)[6]*100, "percent"))
#Ethnicity
print(paste("White/caucasian:", round(prop.table(table(df$race)), 2)[1]*100, "percent"))

##########################################
########################################## 
## Factor Analysis, etc.: Need for Chaos 
## (See Supplementary Materials SM2)
##########################################
########################################## 
nfc <- df %>% 
  select(needchaos_1:needchaos_10) %>% 
  psych::alpha(.) 
print(paste("Cronbach's alpha, Need for Chaos:" , round(nfc$total$raw_alpha,2)))  

#Confirmatory factor analysis, Need for Chaos
cfa <- df %>% 
  select(needchaos_1:needchaos_10) %>% 
  na.omit()

#One-factor (Table S7)
m1  <- ' nfc_oneDim  =~ needchaos_1 + needchaos_3 + needchaos_4 + needchaos_5 + needchaos_6 + needchaos_9 + needchaos_10 
        needchaos_4~~needchaos_5
        needchaos_4~~needchaos_6
        needchaos_5~~needchaos_6'
m1 <- lavaan::cfa(m1, data=cfa) 
m1_fit <- summary(m1 , fit.measures=TRUE) 

print(paste("One-factor model, Chi-square:" , round(m1_fit$FIT["chisq"],2 )))
print(paste("One-factor model, CFI:" , round(m1_fit$FIT["cfi"],2 )))
print(paste("One-factor model, TLI:" , round(m1_fit$FIT["tli"],2 )))
print(paste("One-factor model, RMSEA:" , round(m1_fit$FIT["rmsea"],2 )))

##########################################
########################################## 
## Remove 5% with highest Need for Chaos? (SM8)
## When set to TRUE, used for producing Fig S8, Fig S9, Fig S10
##########################################
##########################################
# Set to TRUE for excluding top 5% of participants with highest Need for Chaos
exclude_obs <- FALSE

if (exclude_obs == TRUE) {
  
  df <- df %>% 
    filter(chaos_z_score < quantile(chaos_z_score, probs = .95, na.rm = TRUE))
  
  print("Excluding top 5% observations with highest levels of the Need for Chaos")
  
}

##########################################
########################################## 
## Analysis - Fig. 3 of Main Text
##########################################
########################################## 
m_share_dem <- lm(share_DEMr_01*-1 ~ chaos_z_score + repdem + female + edu_z_score + age_z_score + race , data = df )
share_dem <- tidy(m_share_dem) %>% 
  filter(term == "chaos_z_score" | term == "repdem") %>% 
  select("term","estimate", "std.error") %>% 
  mutate(dv = "share_dem",
         study = "study4")

m_share_rep <- lm(share_REPr_01*-1 ~ chaos_z_score + repdem + female + edu_z_score + age_z_score + race , data = df )
share_rep <- tidy(m_share_rep) %>% 
  filter(term == "chaos_z_score" | term == "repdem") %>% 
  select("term","estimate", "std.error") %>% 
  mutate(dv = "share_rep",
         study = "study4")

m_true_dem <- lm(true_DEMr_01*-1 ~ chaos_z_score + repdem + female + edu_z_score + age_z_score + race , data = df )
true_dem <- tidy(m_true_dem) %>% 
  filter(term == "chaos_z_score" | term == "repdem") %>% 
  select("term","estimate", "std.error") %>% 
  mutate(dv = "true_dem",
         study = "study4")

m_true_rep <- lm(true_REPr_01*-1 ~ chaos_z_score + repdem + female + edu_z_score + age_z_score + race , data = df )
true_rep <- tidy(m_true_rep) %>% 
  filter(term == "chaos_z_score" | term == "repdem") %>% 
  select("term","estimate", "std.error") %>% 
  mutate(dv = "true_rep",
         study = "study4")

study4_coef <- rbind(share_dem, share_rep, true_dem, true_rep)

#save coefficients for figure 
save(study4_coef,file=here("Data_recoded","study4_fig3.R"))

#regression table for supplemental materials (Table S11)
stargazer(m_share_dem , m_share_rep , m_true_dem , m_true_rep , 
          align=TRUE,
          no.space=TRUE,
          column.labels=c("Share - Dem", "Share - Rep" , "True - Dem" , "True - Rep"),
          covariate.labels = c("Need for Chaos (std.)" , "Party ID (1 = Dem.)" , "Women" , "Education" , "Age  (std.)" , "Ethnicity (1 = non-white)" ,  "Intercept"),
          keep.stat = c("n","ser"),
          dep.var.labels.include = FALSE)

##########################################
########################################## 
## Analysis - Fig. 4 of Main Text
##########################################
########################################## 
#Fig 4: Democratic Voters 
df_dem <- df %>% 
  filter(repdem == 1 )

m_share_dem <- lm(share_DEMr_01*-1 ~ chaos_z_score + female + edu_z_score + age_z_score + race , data = df_dem )
share_dem <- tidy(m_share_dem) %>% 
  filter(term == "chaos_z_score" ) %>% 
  select("term","estimate", "std.error") %>% 
  mutate(dv = "share_dem",
         study = "study4")

m_share_rep <- lm(share_REPr_01*-1 ~ chaos_z_score + female +  edu_z_score + age_z_score + race , data = df_dem )
share_rep <- tidy(m_share_rep) %>% 
  filter(term == "chaos_z_score" ) %>% 
  select("term","estimate", "std.error") %>% 
  mutate(dv = "share_rep",
         study = "study4")

m_true_dem <- lm(true_DEMr_01*-1 ~ chaos_z_score +  female +  edu_z_score + age_z_score + race , data = df_dem )
true_dem <- tidy(m_true_dem) %>% 
  filter(term == "chaos_z_score" ) %>% 
  select("term","estimate", "std.error") %>% 
  mutate(dv = "true_dem",
         study = "study4")

m_true_rep <- lm(true_REPr_01*-1 ~ chaos_z_score + female +  edu_z_score + age_z_score + race , data = df_dem )
true_rep <- tidy(m_true_rep) %>% 
  filter(term == "chaos_z_score" ) %>% 
  select("term","estimate", "std.error") %>% 
  mutate(dv = "true_rep",
         study = "study4")

study4_coef_dem <- rbind(share_dem, share_rep, true_dem, true_rep)

#Fig 4: Republican Voters
df_rep <- df %>% 
  filter(repdem == 0 )

m_share_dem <- lm(share_DEMr_01*-1 ~ chaos_z_score +  female + edu_z_score + age_z_score + race , data = df_rep )
share_dem <- tidy(m_share_dem) %>% 
  filter(term == "chaos_z_score" ) %>% 
  select("term","estimate", "std.error") %>% 
  mutate(dv = "share_dem",
         study = "study4")

m_share_rep <- lm(share_REPr_01*-1 ~ chaos_z_score +  female +  edu_z_score + age_z_score + race , data = df_rep )
share_rep <- tidy(m_share_rep) %>% 
  filter(term == "chaos_z_score" ) %>% 
  select("term","estimate", "std.error") %>% 
  mutate(dv = "share_rep",
         study = "study4")

m_true_dem <- lm(true_DEMr_01*-1 ~ chaos_z_score + female +  edu_z_score + age_z_score + race , data = df_rep )
true_dem <- tidy(m_true_dem) %>% 
  filter(term == "chaos_z_score" ) %>% 
  select("term","estimate", "std.error") %>% 
  mutate(dv = "true_dem",
         study = "study4")

m_true_rep <- lm(true_REPr_01*-1 ~ chaos_z_score +  female +  edu_z_score + age_z_score + race , data = df_rep )
true_rep <- tidy(m_true_rep) %>% 
  filter(term == "chaos_z_score" ) %>% 
  select("term","estimate", "std.error") %>% 
  mutate(dv = "true_rep",
         study = "study4")

study4_coef_rep <- rbind(share_dem, share_rep, true_dem, true_rep)

#save coefficients for figure 
save(study4_coef_dem,file=here("Data_recoded","study4_fig4_dem.R"))
save(study4_coef_rep,file=here("Data_recoded","study4_fig4_rep.R"))

##########################################
########################################## 
## Analysis - Fig. S2 of Supplementary Materials
##########################################
########################################## 
df_independents <- df %>% 
  filter(independents == 1 )

m_share_dem <- lm(share_DEMr_01*-1 ~ chaos_z_score +  female + edu_z_score + age_z_score + race , data = df_independents )
share_dem <- tidy(m_share_dem) %>% 
  filter(term == "chaos_z_score" ) %>% 
  select("term","estimate", "std.error") %>% 
  mutate(dv = "share_dem",
         study = "study4")

m_share_rep <- lm(share_REPr_01*-1 ~ chaos_z_score +   female +  edu_z_score + age_z_score + race , data = df_independents )
share_rep <- tidy(m_share_rep) %>% 
  filter(term == "chaos_z_score" ) %>% 
  select("term","estimate", "std.error") %>% 
  mutate(dv = "share_rep",
         study = "study4")

m_true_dem <- lm(true_DEMr_01*-1 ~ chaos_z_score +   female +  edu_z_score + age_z_score + race , data = df_independents )
true_dem <- tidy(m_true_dem) %>% 
  filter(term == "chaos_z_score" ) %>% 
  select("term","estimate", "std.error") %>% 
  mutate(dv = "true_dem",
         study = "study4")

m_true_rep <- lm(true_REPr_01*-1 ~ chaos_z_score +   female +  edu_z_score + age_z_score + race , data = df_independents )
true_rep <- tidy(m_true_rep) %>% 
  filter(term == "chaos_z_score" ) %>% 
  select("term","estimate", "std.error") %>% 
  mutate(dv = "true_rep",
         study = "study4")

study4_coef_independents <- rbind(share_dem, share_rep, true_dem, true_rep)

#save coefficients for figure 
save(study4_coef_independents,file=here("Data_recoded","study4_figS2.R"))

##########################################
########################################## 
## Analysis - Supplementary Materials SM8
##########################################
########################################## 
#Does Need for Chaos predict intentions to engage in political violence? 
res1 <-  lm(violentact_01 ~ chaos_z_score + repdem  +  female + edu_z_score + age_z_score + race , data = df )

#Table S29
stargazer(res1 , 
          align=TRUE,
          no.space=TRUE,
          column.labels=c("Political Violence Intentions"),
          covariate.labels = c("Need for Chaos (std.)" , "Party ID (1 = Dem.)" , "Women" , "Education" , "Age  (std.)" , "Ethnicity (1 = non-white)" ,  "Intercept"),
          keep.stat = c("n","ser"),
          dep.var.labels.include = FALSE)

##########################################
########################################## 
## Analysis - Fig. 5 of Main Text
##########################################
########################################## 
#Mobilization
mobilizeAll_01 <- lm(mobilizeAll_01 ~ chaos_z_score + repdem  +  female + edu_z_score + age_z_score + race , data = df )
#Amuse
amuseAll_01 <- lm(amuseAll_01 ~ chaos_z_score + repdem  +  female +  edu_z_score + age_z_score + race , data = df )
#Consequences
conseqAll_01 <- lm(conseqAll_01 ~ chaos_z_score + repdem  +  female +  edu_z_score + age_z_score + race , data = df )
#True
trueAll_01 <- lm(trueAll_01 ~ chaos_z_score + repdem  +  female +  edu_z_score + age_z_score + race , data = df )

#Construct data set for plotting
mobilize <- tidy(mobilizeAll_01)[2,]
mobilize$`Motivation for sharing` <- "Mobilize"
mobilize[,1] <- "nfc_mob"

amuse <- tidy(mobilizeAll_01)[2,]
amuse$`Motivation for sharing` <- "Amuse"
amuse[,1] <- "nfc_amuse"

true <- tidy(trueAll_01)[2,]
true$`Motivation for sharing` <- "True"
true[,1] <- "nfc_true"

conse <- tidy(conseqAll_01)[2,]
conse$`Motivation for sharing` <- "Consequences"
conse[,1] <- "nfc_conse"

df_plot <- rbind(mobilize , amuse , true , conse)

text <- paste("beta == ", round(df_plot$estimate , digits = 2) )
df_plot$model <- factor(df_plot$`Motivation for sharing` , levels = c("True" , "Mobilize" , "Consequences" ,  "Amuse"   ))

#Create Figure 5
fig5 <- ggplot(df_plot) +
  geom_hline(yintercept=0 , color = "black") +
  geom_point(aes(x = term , y = estimate, shape = `Motivation for sharing`), size = 3) +
  geom_linerange(aes(x = term, y = estimate, ymin = estimate-2*std.error, ymax = estimate+2*std.error), 
                 size = 2, alpha = .2) +
  theme_bw() +
  ggtitle("Need for Chaos and Motivations for Rumor Sharing" ) +
  ylab("Regression Coefficient Estimate") +   
  xlab("") +
  scale_color_manual(values = c("azure3", "grey45", "grey20" , "black")) +  
  theme(
    axis.ticks.y = element_blank(),
    axis.text.y = element_blank(),
    legend.text=element_text(size=13),
    axis.text.x = element_text(size = 13 , color = "black"),
    axis.title.x = element_text(size = 14)) + 
  geom_text(aes(x = term , y = estimate),
    label=text, 
    nudge_y = 0.0, nudge_x = .2, 
    check_overlap = T , parse = TRUE , size = 6
  ) + 
  coord_flip() + theme(legend.justification=c(1,0), legend.position=c(.99,.67)) +
  labs(color = "Motivation for sharing") +
  guides(fill = guide_legend(reverse=TRUE), shape = guide_legend(reverse = TRUE)) +
  scale_y_continuous(breaks = c(-.05,-.025,0,.025,.05,.075,.1), 
                     labels = c("-.05", "-.025", "0", ".025", ".05", ".075" ,".1"),
                     limits = c(-.05,.1))

fig5

ggsave(here::here("figures", "fig5.pdf"),
       width = 8, height = 5)

#regression table for Supplemental Materials (Table S15)
stargazer(mobilizeAll_01 , amuseAll_01 , conseqAll_01 , trueAll_01 , 
          align=TRUE,
          no.space=TRUE,
          column.labels=c("Mobilize", "Amuse" , "Conseq." , "True"),
          covariate.labels = c("Need for Chaos" , "Party Id (1 = Dem)" , "Gender (1 = female)" ,"Education","Age", "Ethnicity (1 = non-white) ", "Intercept"),
          keep.stat = c("n","ser"),
          dep.var.labels.include = FALSE
)

#regression table for Supplemental Materials (Table S16)
mobilizeDem_01 <- lm(mobilizeDem_01 ~ chaos_z_score + repdem  +  female + edu_z_score + age_z_score + race , data = df )
mobilizeRep_01 <- lm(mobilizeRep_01 ~ chaos_z_score + repdem  +  female +  edu_z_score + age_z_score + race , data = df )

amuseDem_01 <- lm(amuseDem_01 ~ chaos_z_score + repdem  +  female +  edu_z_score + age_z_score + race , data = df )
amuseRep_01 <- lm(amuseRep_01 ~ chaos_z_score + repdem  +  female +  edu_z_score + age_z_score + race , data = df )

conseqDem_01 <- lm(conseqDem_01 ~ chaos_z_score + repdem  +  female +  edu_z_score + age_z_score + race , data = df )
conseqRep_01 <- lm(conseqRep_01 ~ chaos_z_score + repdem  +  female +  edu_z_score + age_z_score + race , data = df )

trueDem_01 <- lm(trueDem_01 ~ chaos_z_score + repdem  +  female +  edu_z_score + age_z_score + race , data = df )
trueRep_01 <- lm(trueRep_01 ~ chaos_z_score + repdem  +  female +  edu_z_score + age_z_score + race , data = df )

stargazer(mobilizeDem_01 , mobilizeRep_01 , amuseDem_01 , amuseRep_01 , conseqDem_01, conseqRep_01, trueDem_01, trueRep_01,
          align=TRUE,
          no.space=TRUE,
          column.labels=c("Mob. Dem", "Mob. Rep" , "Amuse Dem" , "Amuse Rep", "Consq. Dem", "Consq. Rep", "True Dem", "True Rep"),
          covariate.labels = c("Need for Chaos" , "Party Id (1 = Dem)" , "Women" ,"Education","Age", "Ethnicity (1 = non-white) ", "Intercept"),
          keep.stat = c("n"),
          dep.var.labels.include = FALSE
)
